rm(list=ls())

library("ggplot2")
library("foreign")
library('ggpubr')
library('gridExtra')

ExpR <- read.csv("matlab/baseline_model/conditionalInterventionProbs_replication.csv")
Lake <- read.dta("r/Lake_HIR_Country_yearreplication.dta")
Lake <- subset(Lake, year <= 1999)
LakeSH <- aggregate(Lake$us_SH1995, list(Lake$ccode), mean, na.rm=T)

ExpR <- merge(ExpR, LakeSH, by.x="ccode", by.y = "Group.1")

ExpR$Warsaw <- rep(0,dim(ExpR)[1])
ExpR$Warsaw[ExpR$ccode %in% c(290, 300, 310, 315, 316, 317, 
                              339, 355, 359, 360, 
                              366, 367, 368, 369, 370, 371, 372, 373, 
                              640, 701, 702, 703, 704, 705)] <- 1

### Warsaw expectations
# mean(ExpR$RUS[ExpR$Warsaw==1])
# mean(ExpR$RUS[ExpR$Warsaw==0])
# t.test(RUS ~ Warsaw, data=ExpR)

# highest warsaw is HUN
# highest no warsaw is FIN, Qatar, UAE, and North Korea
xalt <- 0.05
ExpR$WarsawW <- factor(ifelse(ExpR$Warsaw==1,"Member", "Nonmember"), levels=c("Nonmember","Member"))
pbox <- ggplot(ExpR, aes(y=RUS, x=WarsawW)) + geom_boxplot() +  
  theme_classic(18) + ylab("RUS Intervention Probability") + xlab("Warsaw Pact Membership") +
  theme(axis.title.x=element_text(margin = margin(t = 15, r = 0, b = 0, l = 0))) +
  theme(axis.title.y=element_text(margin = margin(t = 0, r = 15, b = 0, l = 0))) +
  annotate("text", x = 2.05+xalt, y = 0.585, label = "HUN") +
  annotate("text", x = 1.05+xalt, y = 0.624, label = "QAT") + 
  annotate("text", x = 1.05+xalt, y = 0.567, label = "FIN") + 
  annotate("text", x = 0.95-xalt, y = 0.562, label = "UAE") + 
  annotate("text", x = 1.05+xalt, y = 0.555, label = "PRK") +
  scale_y_continuous(breaks = seq(from=0.4, to=0.65, by=0.05), labels=100*seq(from=0.4, to=0.65, by=0.05)) 


xalt <- 0.025
psct <- ggplot(ExpR, aes(x=x, y=US)) + geom_point(size=1.5)+ geom_smooth(method=lm, color="black") +
  theme_classic(18)+ xlab("Lake's (2009) U.S. Security Hierarchy") + ylab("U.S. Intervention Probability") +
  theme(axis.title.x=element_text(margin = margin(t = 15, r = 0, b = 0, l = 0))) +
  theme(axis.title.y=element_text(margin = margin(t = 0, r = 15, b = 0, l = 0))) +
  annotate("text", x = 0.05+xalt, y = 0.6, label = "PRK") + 
  annotate("text", x = 0.05+xalt, y = 0.535, label = "IRQ") +
  annotate("text", x = 0.05+xalt, y = 0.57, label = "QAT") +
  annotate("text", x = 1.57-xalt, y = 0.531, label = "PAN") + 
  annotate("text", x = 1.19+xalt, y = 0.47, label = "VNM") +
  annotate("text", x = 0.05+xalt, y = 0.394, label = "LKA") +
  annotate("text", x = 0.05+xalt, y = 0.401, label = "BGD") + 
  annotate("text", x = 0.6+xalt, y = 0.5825, label = "CAN") + 
  annotate("text", x = 0.665+xalt, y = 0.47, label = "JPN") + 
  annotate("text", x = 0.415+xalt, y = 0.44, label = "NZL") + 
  annotate("text", x = 0.55+xalt, y = 0.5275, label = "HND")+
  scale_y_continuous(breaks = seq(from=0.4, to=0.65, by=0.05), labels=100*seq(from=0.4, to=0.65, by=0.05))

finalgraph <- arrangeGrob(pbox, psct, nrow = 1)
ggsave('fig1.pdf', plot = finalgraph,width =14, height=7)